import numpy as np
import matplotlib.pyplot as plt

plt.rcParams['font.sans-serif'] = ['SimHei']  # 用黑体显示中文
plt.rcParams['axes.unicode_minus'] = False    # 正常显示负号

# 实验三：复变函数的级数展开 (简化版)
print("实验三：复变函数的级数展开")
print("="*40)

# 计算幂级数的部分和
def power_series_sum(coefficients, z, n_terms=50):
    z_power = 1
    series_sum = 0
    for i in range(min(n_terms, len(coefficients))):
        series_sum += coefficients[i] * z_power
        z_power *= z
    return series_sum

# (1-z)^(-m) 的幂级数展开（习题21）
def binomial_series_coeffs(m, n_terms=20):
    coeffs = []
    for n in range(n_terms):
        coeff = np.math.factorial(n + m - 1) // (np.math.factorial(m - 1) * np.math.factorial(n))
        coeffs.append(coeff)
    return coeffs

print("1. (1-z)^(-m) 的幂级数展开:")
m = 3
coeffs = binomial_series_coeffs(m, 10)
print(f"  (1-z)^(-{m}) 的前10项系数:")
for n, coeff in enumerate(coeffs):
    print(f"    a_{n} = {coeff}")

# 验证在收敛圆盘内级数的收敛性
z_test = 0.5  # |z| < 1，应在收敛圆盘内
approximation = power_series_sum(coeffs, z_test, 20)
exact_value = (1 - z_test) ** (-m)
print(f"\n在 z = {z_test} 处:")
print(f"  级数近似值: {approximation:.6f}")
print(f"  精确值: {exact_value:.6f}")
print(f"  误差: {abs(approximation - exact_value):.2e}")

# (2z+3)/(z+1) 展开为 (z-1) 的幂级数（习题22）
print(f"\n2. (2z+3)/(z+1) 展开为 (z-1) 的幂级数:")
def rational_series_coeffs(n_terms=10):
    coeffs = [2.5]  # 第一项系数
    for n in range(1, n_terms):
        coeffs.append((-1)**n / (2**(n+1)))
    return coeffs

coeffs2 = rational_series_coeffs(8)
print(f"  前8项系数:")
for n, coeff in enumerate(coeffs2):
    print(f"    a_{n} = {coeff:.6f}")
print(f"  收敛半径: 2")

# 验证级数
z_test2 = 1.5  # z-1 = 0.5, |z-1| = 0.5 < 2
w_test = z_test2 - 1  # w = z-1
approximation2 = power_series_sum(coeffs2, w_test, 20)
exact_value2 = (2*z_test2 + 3) / (z_test2 + 1)
print(f"\n在 z = {z_test2} 处:")
print(f"  级数近似值: {approximation2:.6f}")
print(f"  精确值: {exact_value2:.6f}")
print(f"  误差: {abs(approximation2 - exact_value2):.2e}")

# 可视化级数收敛性
print("\n3. 可视化 (1-z)^(-2) 的级数部分和:")
z_real = np.linspace(-0.8, 0.8, 80)
partial_sums_10 = [power_series_sum(binomial_series_coeffs(2, 20), z) for z in z_real]
partial_sums_20 = [power_series_sum(binomial_series_coeffs(2, 30), z) for z in z_real]
exact_vals = [(1-z)**(-2) for z in z_real]

fig, ax = plt.subplots(1, 1, figsize=(8, 5))
ax.plot(z_real, exact_vals, 'k-', label='精确值', linewidth=2)
ax.plot(z_real, [np.real(s) for s in partial_sums_10], 'r--', label='前10项', linewidth=1)
ax.plot(z_real, [np.real(s) for s in partial_sums_20], 'b-.', label='前20项', linewidth=1)
ax.set_title('$(1-z)^{(-2)}$ 的幂级数展开收敛性')
ax.grid(True)
ax.legend()
plt.tight_layout()
plt.show()